#######################################################
#####REPLICATION ARCHIVE FOR NUANCED ACCOUNTABILITY####
#########BRITISH JOURNAL OF POLITICAL SCIENCE##########
###########WARD DATA, LOCAL ELECTIONS ONLY#############
###de Kadt & Lieberman, 2017

#######################################################
################ SETUP BEGINS HERE ####################
#######################################################
## USER: change working directory on line 56 ##########
###Packages
library(foreign)
library(stargazer)
library(Hmisc)
library(sandwich)
library(lmtest)

###Cluster Robust SEs###
vcovCluster <- function(
  model,
  cluster
)
{
  require(sandwich)
  require(lmtest)
  if(nrow(model.matrix(model))!=length(cluster)){
    stop("check your data: cluster variable has different N than model")
  }
  M <- length(unique(cluster))
  N <- length(cluster)           
  K <- model$rank   
  if(M<50){
    warning("Fewer than 50 clusters, variances may be unreliable (could try block bootstrap instead).")
  }
  dfc <- (M/(M - 1)) * ((N - 1)/(N - K))
  uj  <- apply(estfun(model), 2, function(x) tapply(x, cluster, sum));
  rcse.cov <- dfc * sandwich(model, meat = crossprod(uj)/N)
  return(rcse.cov)
}

###Correlation Table###
cor.prob <- function(X, dfr = nrow(X) - 2) {
  R <- cor(X)
  above <- row(R) < col(R)
  r2 <- R[above]^2
  Fstat <- r2 * dfr / (1 - r2)
  R[above] <- 1 - pf(Fstat, 1, dfr)
  
  cor.mat <- t(R)
  cor.mat[upper.tri(cor.mat)] <- NA
  cor.mat
}

###Working Directories###
## USER: Input your own working directory for files on the next line ##
dir = paste("####INPUT WORKING FILE DIRECTORY BETWEEN QUOTATION MARKS ####")

##

fig = paste(dir,"output/", sep="")
stgz =  paste(dir,"output/", sep="")

###Data Setup###
#Load dataset
dat = read.dta(paste(dir,"nuanced_accountability_ward.dta", sep=""))
#Set up covariates
cov = paste("+ unempfrac_ch10 + log_income_ch10 + log_pop_ch10 + femalefrac_ch10 + whitefrac_ch10 + employment_broad2001 + log_income2001 + log_pop2001 + sex2001  + white_frac2001 + share_area_all_tbvc_ta + opp_control")

#######################################################
################ MAIN RESULTS HERE ####################
#######################################################

###################################
####LOCAL GOVERNMENT ELECTIONS#####
###################################
i = 1
#General holders for results:
out.mod.full.nocov = list()
out.mod.full.cov = list()
out.mod.full.nocov.comp = list()
out.mod.full.cov.comp = list()
out.mod.full.interact.cov = list()
out.mod.full.interact.cov.comp = list()

out.mod.full.nocov.incumb = list()
out.mod.full.nocov.incumb.comp = list()
out.mod.full.cov.incumb = list()
out.mod.full.cov.incumb.comp = list()
out.mod.full.nocov.opp = list()
out.mod.full.nocov.opp.comp = list()
out.mod.full.cov.opp = list()
out.mod.full.cov.opp.comp = list()

out.mod.full.cov.nowc = list()
out.mod.full.cov.comp.nowc = list()
out.mod.full.cov.incumb.nowc = list()
out.mod.full.cov.incumb.comp.nowc = list()
out.mod.full.cov.opp.nowc = list()
out.mod.full.cov.opp.comp.nowc = list()

out.mod.full.cov.nokzn = list()
out.mod.full.cov.comp.nokzn = list()
out.mod.full.cov.incumb.nokzn = list()
out.mod.full.cov.incumb.comp.nokzn = list()
out.mod.full.cov.opp.nokzn = list()
out.mod.full.cov.opp.comp.nokzn = list()

out.mod.full.cov.noec = list()
out.mod.full.cov.comp.noec = list()
out.mod.full.cov.incumb.noec = list()
out.mod.full.cov.incumb.comp.noec = list()
out.mod.full.cov.opp.noec = list()
out.mod.full.cov.opp.comp.noec = list()

lm.out.full.nocov = list()
lm.out.full.cov = list()
lm.out.full.nocov.comp = list()
lm.out.full.cov.comp = list()
lm.out.full.interact.cov = list()
lm.out.full.interact.cov.comp = list()

lm.out.full.nocov.incumb = list()
lm.out.full.nocov.incumb.comp = list()
lm.out.full.cov.incumb = list()
lm.out.full.cov.incumb.comp = list()
lm.out.full.nocov.opp = list()
lm.out.full.nocov.opp.comp = list()
lm.out.full.cov.opp = list()
lm.out.full.cov.opp.comp = list()

lm.out.full.cov.nowc = list()
lm.out.full.cov.comp.nowc = list()
lm.out.full.cov.incumb.nowc = list()
lm.out.full.cov.incumb.comp.nowc = list()
lm.out.full.cov.opp.nowc = list()
lm.out.full.cov.opp.comp.nowc = list()

lm.out.full.cov.nokzn = list()
lm.out.full.cov.comp.nokzn = list()
lm.out.full.cov.incumb.nokzn = list()
lm.out.full.cov.incumb.comp.nokzn = list()
lm.out.full.cov.opp.nokzn = list()
lm.out.full.cov.opp.comp.nokzn = list()

lm.out.full.cov.noec = list()
lm.out.full.cov.comp.noec = list()
lm.out.full.cov.incumb.noec = list()
lm.out.full.cov.incumb.comp.noec = list()
lm.out.full.cov.opp.noec = list()
lm.out.full.cov.opp.comp.noec = list()

for(iv in c("water_ch10", "flush_ch10", "refuse_ch10", "sd_means_ch10")){
  dat.80th = dat[which(dat$water_80th==1 & dat$flush_80th==1 & dat$refuse_80th==1),]
  
  dat.sub = dat.80th[,c(iv,
                        "lgeanc_vs_ch10", "unempfrac_ch10", "log_income_ch10", "log_pop_ch10", "femalefrac_ch10", "blackfrac_ch10", 
                        "whitefrac_ch10", "indianfrac_ch10", "colouredfrac_ch10", "traddwell_ch10", "cat_b", "opp_control", "province2011",
                        "employment_broad2001","log_income2001","log_pop2001","sex2001","indian_frac2001","colored_frac2001","black_frac2001","white_frac2001",
                        "share_area_all_tbvc_ta","traddwell_proportion2001","refuse_perc2001","all_flush_perc2001","home_water_perc2001","sd_means2001",
                        "competitive")]
  dat.sub = na.omit(dat.sub)
  
  dat.sub.comp = dat.sub[which(dat.sub$competitive==1),]
  
  dat.sub.comp = na.omit(dat.sub.comp)
  
  iv.int = paste(iv,"*opp_control", sep="")
  
  if(iv=="water_ch10"){
    baseline = " + home_water_perc2001"
  } else if(iv=="flush_ch10"){
    baseline = " + all_flush_perc2001"
  } else if(iv=="refuse_ch10"){
    baseline = " + refuse_perc2001"
  } else {
    baseline = " + sd_means2001"
  }  
  
  #####REGRESSIONS:
  ##Main specifications
  lm.out.full.nocov[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv, sep=""), data=dat.sub)
  lm.out.full.cov[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv, cov, baseline, sep=""), data=dat.sub)
  
  lm.out.full.nocov.comp[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv, sep=""), data=dat.sub.comp)
  lm.out.full.cov.comp[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv, cov, baseline, sep=""), data=dat.sub.comp)
  
  lm.out.full.interact.cov[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub)
  lm.out.full.interact.cov.comp[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub.comp)
  
  ##By incumbency status
  lm.out.full.nocov.incumb[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, sep=""), data=dat.sub[which(dat.sub$opp_control==0),])
  lm.out.full.cov.incumb[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub[which(dat.sub$opp_control==0),])
  lm.out.full.nocov.incumb.comp[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, sep=""), data=dat.sub.comp[which(dat.sub.comp$opp_control==0),])
  lm.out.full.cov.incumb.comp[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub.comp[which(dat.sub.comp$opp_control==0),])
  
  lm.out.full.nocov.opp[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, sep=""), data=dat.sub[which(dat.sub$opp_control==1),])
  lm.out.full.cov.opp[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub[which(dat.sub$opp_control==1),])
  lm.out.full.nocov.opp.comp[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, sep=""), data=dat.sub.comp[which(dat.sub.comp$opp_control==1),])
  lm.out.full.cov.opp.comp[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub.comp[which(dat.sub.comp$opp_control==1),])
  
  ##Excluding particular provinces
  lm.out.full.cov.nowc[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub[which(dat.sub$province2011!="Western Cape"),])
  lm.out.full.cov.comp.nowc[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub.comp[which(dat.sub.comp$province2011!="Western Cape"),])
  lm.out.full.cov.incumb.nowc[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub[which(dat.sub$opp_control==0 & dat.sub$province2011!="Western Cape"),])
  lm.out.full.cov.incumb.comp.nowc[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub.comp[which(dat.sub.comp$opp_control==0 & dat.sub.comp$province2011!="Western Cape"),])
  lm.out.full.cov.opp.nowc[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub[which(dat.sub$opp_control==1 & dat.sub$province2011!="Western Cape"),])
  lm.out.full.cov.opp.comp.nowc[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub.comp[which(dat.sub.comp$opp_control==1 & dat.sub.comp$province2011!="Western Cape"),])
  
  lm.out.full.cov.nokzn[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub[which(dat.sub$province2011!="KwaZulu-Natal"),])
  lm.out.full.cov.comp.nokzn[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub.comp[which(dat.sub.comp$province2011!="KwaZulu-Natal"),])
  lm.out.full.cov.incumb.nokzn[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub[which(dat.sub$opp_control==0 & dat.sub$province2011!="KwaZulu-Natal"),])
  lm.out.full.cov.incumb.comp.nokzn[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub.comp[which(dat.sub.comp$opp_control==0 & dat.sub.comp$province2011!="KwaZulu-Natal"),])
  lm.out.full.cov.opp.nokzn[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub[which(dat.sub$opp_control==1 & dat.sub$province2011!="KwaZulu-Natal"),])
  lm.out.full.cov.opp.comp.nokzn[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub.comp[which(dat.sub.comp$opp_control==1 & dat.sub.comp$province2011!="KwaZulu-Natal"),])
  
  lm.out.full.cov.noec[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub[which(dat.sub$province2011!="Eastern Cape"),])
  lm.out.full.cov.comp.noec[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub.comp[which(dat.sub.comp$province2011!="Eastern Cape"),])
  lm.out.full.cov.incumb.noec[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub[which(dat.sub$opp_control==0 & dat.sub$province2011!="Eastern Cape"),])
  lm.out.full.cov.incumb.comp.noec[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub.comp[which(dat.sub.comp$opp_control==0 & dat.sub.comp$province2011!="Eastern Cape"),])
  lm.out.full.cov.opp.noec[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub[which(dat.sub$opp_control==1 & dat.sub$province2011!="Eastern Cape"),])
  lm.out.full.cov.opp.comp.noec[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub.comp[which(dat.sub.comp$opp_control==1 & dat.sub.comp$province2011!="Eastern Cape"),])
  
  #######CLUSTER STANDARD ERRORS:
  out.mod.full.nocov[[i]] = coeftest(lm.out.full.nocov[[i]],vcov = vcovCluster(lm.out.full.nocov[[i]], dat.sub$cat_b))
  out.mod.full.cov[[i]] = coeftest(lm.out.full.cov[[i]],vcov = vcovCluster(lm.out.full.cov[[i]], dat.sub$cat_b))
  
  out.mod.full.nocov.comp[[i]] = coeftest(lm.out.full.nocov.comp[[i]],vcov = vcovCluster(lm.out.full.nocov.comp[[i]], dat.sub.comp$cat_b))
  out.mod.full.cov.comp[[i]] = coeftest(lm.out.full.cov.comp[[i]],vcov = vcovCluster(lm.out.full.cov.comp[[i]], dat.sub.comp$cat_b))
  
  out.mod.full.interact.cov[[i]] = coeftest(lm.out.full.interact.cov[[i]],vcov = vcovCluster(lm.out.full.interact.cov[[i]], dat.sub$cat_b))
  out.mod.full.interact.cov.comp[[i]] = coeftest(lm.out.full.interact.cov.comp[[i]],vcov = vcovCluster(lm.out.full.interact.cov.comp[[i]], dat.sub.comp$cat_b))
  
  out.mod.full.nocov.incumb[[i]] = coeftest(lm.out.full.nocov.incumb[[i]],vcov = vcovCluster(lm.out.full.nocov.incumb[[i]], dat.sub[which(dat.sub$opp_control==0),]$cat_b))
  out.mod.full.cov.incumb[[i]] = coeftest(lm.out.full.cov.incumb[[i]],vcov = vcovCluster(lm.out.full.cov.incumb[[i]], dat.sub[which(dat.sub$opp_control==0),]$cat_b))
  out.mod.full.nocov.incumb.comp[[i]] = coeftest(lm.out.full.nocov.incumb.comp[[i]],vcov = vcovCluster(lm.out.full.nocov.incumb.comp[[i]], dat.sub.comp[which(dat.sub.comp$opp_control==0),]$cat_b))
  out.mod.full.cov.incumb.comp[[i]] = coeftest(lm.out.full.cov.incumb.comp[[i]],vcov = vcovCluster(lm.out.full.cov.incumb.comp[[i]], dat.sub.comp[which(dat.sub.comp$opp_control==0),]$cat_b))
  
  out.mod.full.nocov.opp[[i]] = coeftest(lm.out.full.nocov.opp[[i]],vcov = vcovCluster(lm.out.full.nocov.opp[[i]], dat.sub[which(dat.sub$opp_control==1),]$cat_b))
  out.mod.full.cov.opp[[i]] = coeftest(lm.out.full.cov.opp[[i]],vcov = vcovCluster(lm.out.full.cov.opp[[i]], dat.sub[which(dat.sub$opp_control==1),]$cat_b))
  out.mod.full.nocov.opp.comp[[i]] = coeftest(lm.out.full.nocov.opp.comp[[i]],vcov = vcovCluster(lm.out.full.nocov.opp.comp[[i]], dat.sub.comp[which(dat.sub.comp$opp_control==1),]$cat_b))
  out.mod.full.cov.opp.comp[[i]] = coeftest(lm.out.full.cov.opp.comp[[i]],vcov = vcovCluster(lm.out.full.cov.opp.comp[[i]], dat.sub.comp[which(dat.sub.comp$opp_control==1),]$cat_b))
  
  out.mod.full.cov.nowc[[i]] = coeftest(lm.out.full.cov.nowc[[i]],vcov = vcovCluster(lm.out.full.cov.nowc[[i]], dat.sub[which(dat.sub$province2011!="Western Cape"),]$cat_b))
  out.mod.full.cov.comp.nowc[[i]] = coeftest(lm.out.full.cov.comp.nowc[[i]],vcov = vcovCluster(lm.out.full.cov.comp.nowc[[i]], dat.sub.comp[which(dat.sub.comp$province2011!="Western Cape"),]$cat_b))
  out.mod.full.cov.incumb.nowc[[i]] = coeftest(lm.out.full.cov.incumb.nowc[[i]],vcov = vcovCluster(lm.out.full.cov.incumb.nowc[[i]], dat.sub[which(dat.sub$opp_control==0 & dat.sub$province2011!="Western Cape"),]$cat_b))
  out.mod.full.cov.incumb.comp.nowc[[i]] = coeftest(lm.out.full.cov.incumb.comp.nowc[[i]],vcov = vcovCluster(lm.out.full.cov.incumb.comp.nowc[[i]], dat.sub.comp[which(dat.sub.comp$opp_control==0 & dat.sub.comp$province2011!="Western Cape"),]$cat_b))
  out.mod.full.cov.opp.nowc[[i]] = coeftest(lm.out.full.cov.opp.nowc[[i]],vcov = vcovCluster(lm.out.full.cov.opp.nowc[[i]], dat.sub[which(dat.sub$opp_control==1 & dat.sub$province2011!="Western Cape"),]$cat_b))
  out.mod.full.cov.opp.comp.nowc[[i]] = coeftest(lm.out.full.cov.opp.comp.nowc[[i]],vcov = vcovCluster(lm.out.full.cov.opp.comp.nowc[[i]], dat.sub.comp[which(dat.sub.comp$opp_control==1 & dat.sub.comp$province2011!="Western Cape"),]$cat_b))
  
  out.mod.full.cov.nokzn[[i]] = coeftest(lm.out.full.cov.nokzn[[i]],vcov = vcovCluster(lm.out.full.cov.nokzn[[i]], dat.sub[which(dat.sub$province2011!="KwaZulu-Natal"),]$cat_b))
  out.mod.full.cov.comp.nokzn[[i]] = coeftest(lm.out.full.cov.comp.nokzn[[i]],vcov = vcovCluster(lm.out.full.cov.comp.nokzn[[i]], dat.sub.comp[which(dat.sub.comp$province2011!="KwaZulu-Natal"),]$cat_b))
  out.mod.full.cov.incumb.nokzn[[i]] = coeftest(lm.out.full.cov.incumb.nokzn[[i]],vcov = vcovCluster(lm.out.full.cov.incumb.nokzn[[i]], dat.sub[which(dat.sub$opp_control==0 & dat.sub$province2011!="KwaZulu-Natal"),]$cat_b))
  out.mod.full.cov.incumb.comp.nokzn[[i]] = coeftest(lm.out.full.cov.incumb.comp.nokzn[[i]],vcov = vcovCluster(lm.out.full.cov.incumb.comp.nokzn[[i]], dat.sub.comp[which(dat.sub.comp$opp_control==0 & dat.sub.comp$province2011!="KwaZulu-Natal"),]$cat_b))
  out.mod.full.cov.opp.nokzn[[i]] = coeftest(lm.out.full.cov.opp.nokzn[[i]],vcov = vcovCluster(lm.out.full.cov.opp.nokzn[[i]], dat.sub[which(dat.sub$opp_control==1 & dat.sub$province2011!="KwaZulu-Natal"),]$cat_b))
  out.mod.full.cov.opp.comp.nokzn[[i]] = coeftest(lm.out.full.cov.opp.comp.nokzn[[i]],vcov = vcovCluster(lm.out.full.cov.opp.comp.nokzn[[i]], dat.sub.comp[which(dat.sub.comp$opp_control==1 & dat.sub.comp$province2011!="KwaZulu-Natal"),]$cat_b))
  
  out.mod.full.cov.noec[[i]] = coeftest(lm.out.full.cov.noec[[i]],vcov = vcovCluster(lm.out.full.cov.noec[[i]], dat.sub[which(dat.sub$province2011!="Eastern Cape"),]$cat_b))
  out.mod.full.cov.comp.noec[[i]] = coeftest(lm.out.full.cov.comp.noec[[i]],vcov = vcovCluster(lm.out.full.cov.comp.noec[[i]], dat.sub.comp[which(dat.sub.comp$province2011!="Eastern Cape"),]$cat_b))
  out.mod.full.cov.incumb.noec[[i]] = coeftest(lm.out.full.cov.incumb.noec[[i]],vcov = vcovCluster(lm.out.full.cov.incumb.noec[[i]], dat.sub[which(dat.sub$opp_control==0 & dat.sub$province2011!="Eastern Cape"),]$cat_b))
  out.mod.full.cov.incumb.comp.noec[[i]] = coeftest(lm.out.full.cov.incumb.comp.noec[[i]],vcov = vcovCluster(lm.out.full.cov.incumb.comp.noec[[i]], dat.sub.comp[which(dat.sub.comp$opp_control==0 & dat.sub.comp$province2011!="Eastern Cape"),]$cat_b))
  out.mod.full.cov.opp.noec[[i]] = coeftest(lm.out.full.cov.opp.noec[[i]],vcov = vcovCluster(lm.out.full.cov.opp.noec[[i]], dat.sub[which(dat.sub$opp_control==1 & dat.sub$province2011!="Eastern Cape"),]$cat_b))
  out.mod.full.cov.opp.comp.noec[[i]] = coeftest(lm.out.full.cov.opp.comp.noec[[i]],vcov = vcovCluster(lm.out.full.cov.opp.comp.noec[[i]], dat.sub.comp[which(dat.sub.comp$opp_control==1 & dat.sub.comp$province2011!="Eastern Cape"),]$cat_b))
  
  i = i+1
}

###Generate appropriate SEs for the tables
cl.se.full.nocov = lapply(out.mod.full.nocov, `[`,,2)
cl.se.full.nocov.comp = lapply(out.mod.full.nocov.comp, `[`,,2)
cl.se.full.cov = lapply(out.mod.full.cov, `[`,,2)
cl.se.full.cov.comp = lapply(out.mod.full.cov.comp, `[`,,2)
cl.se.full.interact.cov = lapply(out.mod.full.interact.cov, `[`,,2)
cl.se.full.interact.cov.comp = lapply(out.mod.full.interact.cov.comp, `[`,,2)
cl.se.full.nocov.incumb = lapply(out.mod.full.nocov.incumb, `[`,,2)
cl.se.full.cov.incumb = lapply(out.mod.full.cov.incumb, `[`,,2)
cl.se.full.nocov.incumb.comp = lapply(out.mod.full.nocov.incumb.comp, `[`,,2)
cl.se.full.cov.incumb.comp = lapply(out.mod.full.cov.incumb.comp, `[`,,2)
cl.se.full.cov.opp = lapply(out.mod.full.cov.opp, `[`,,2)
cl.se.full.nocov.opp = lapply(out.mod.full.nocov.opp, `[`,,2)
cl.se.full.cov.opp.comp = lapply(out.mod.full.cov.opp.comp, `[`,,2)
cl.se.full.nocov.opp.comp = lapply(out.mod.full.nocov.opp.comp, `[`,,2)

cl.se.full.cov.comp.nowc = lapply(out.mod.full.cov.comp.nowc, `[`,,2)
cl.se.full.cov.comp.nokzn = lapply(out.mod.full.cov.comp.nokzn, `[`,,2)
cl.se.full.cov.comp.noec = lapply(out.mod.full.cov.comp.noec, `[`,,2)

cl.se.full.cov.nowc = lapply(out.mod.full.cov.nowc, `[`,,2)
cl.se.full.cov.nokzn = lapply(out.mod.full.cov.nokzn, `[`,,2)
cl.se.full.cov.noec = lapply(out.mod.full.cov.noec, `[`,,2)

###Table 4 in the paper (panel combination of the following two stargazer tables): 
#Heterogeneous Effects of Change in SD on ANC Voteshare by ANC incumbency status, LGE
filename = paste("lgecovariatesANC.tex")
stargazer(lm.out.full.cov.incumb[[1]], lm.out.full.cov.incumb.comp[[1]],
          lm.out.full.cov.incumb[[2]], lm.out.full.cov.incumb.comp[[2]],
          lm.out.full.cov.incumb[[3]], lm.out.full.cov.incumb.comp[[3]],
          lm.out.full.cov.incumb[[4]], lm.out.full.cov.incumb.comp[[4]],
          #This next row includes the correct SEs:
          se = list(cl.se.full.cov.incumb[[1]], cl.se.full.cov.incumb.comp[[1]],
                    cl.se.full.cov.incumb[[2]], cl.se.full.cov.incumb.comp[[2]],
                    cl.se.full.cov.incumb[[3]], cl.se.full.cov.incumb.comp[[3]],
                    cl.se.full.cov.incumb[[4]], cl.se.full.cov.incumb.comp[[4]]),
          style = "ajps",
          label = paste("tab:lge-covariates-anconly"),
          column.separate = c(8),
          #covariate.labels = covars,
          title = "Heterogeneous effects of change in service delivery on change on local ANC vote: ANC municipalities only",
          model.numbers=F,
          dep.var.labels.include = FALSE,
          column.labels = c("Change in ANC Vote Share (Local Elections)"),
          #column.separate = c(1,1,1,1,1),
          omit.stat = c("f", "adj.rsq", "ser"),
          notes = "Standard errors clustered by municipality",
          float.env = "sidewaystable",
          font.size = "tiny",
          add.lines = list(c("Covariates?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("$<$80th Percentile SD?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("ANC Control Only?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("Competitive only?","","Y","","Y","","Y","","Y")),
          out = paste(stgz,filename,sep="\\"))

filename = paste("lgecovariatesOPP.tex")
stargazer(lm.out.full.cov.opp[[1]], lm.out.full.cov.opp.comp[[1]],
          lm.out.full.cov.opp[[2]], lm.out.full.cov.opp.comp[[2]],
          lm.out.full.cov.opp[[3]], lm.out.full.cov.opp.comp[[3]],
          lm.out.full.cov.opp[[4]], lm.out.full.cov.opp.comp[[4]],
          #This next row includes the correct SEs:
          se = list(cl.se.full.cov.opp[[1]], cl.se.full.cov.opp.comp[[1]],
                    cl.se.full.cov.opp[[2]], cl.se.full.cov.opp.comp[[2]],
                    cl.se.full.cov.opp[[3]], cl.se.full.cov.opp.comp[[3]],
                    cl.se.full.cov.opp[[4]], cl.se.full.cov.opp.comp[[4]]),
          style = "ajps",
          label = paste("tab:lge-covariates-opponly"),
          column.separate = c(8),
          #covariate.labels = covars,
          title = "Heterogeneous effects of change in service delivery on change on local ANC vote: Opposition municipalities only",
          model.numbers=F,
          dep.var.labels.include = FALSE,
          column.labels = c("Change in ANC Vote Share (Local Elections)"),
          #column.separate = c(1,1,1,1,1),
          omit.stat = c("f", "adj.rsq", "ser"),
          notes = "Standard errors clustered by municipality",
          float.env = "sidewaystable",
          font.size = "tiny",
          add.lines = list(c("Covariates?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("$<$80th Percentile SD?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("Opp Control Only?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("Competitive only?","","Y","","Y","","Y","","Y")),
          out = paste(stgz,filename,sep="\\"))

###APPENDIX TABLES###
###Table 18 in Appendix: LGE, ANC only, No Covariates
filename = paste("lgenocovariatesANC.tex")
stargazer(lm.out.full.nocov.incumb[[1]], lm.out.full.nocov.incumb.comp[[1]],
          lm.out.full.nocov.incumb[[2]], lm.out.full.nocov.incumb.comp[[2]],
          lm.out.full.nocov.incumb[[3]], lm.out.full.nocov.incumb.comp[[3]],
          lm.out.full.nocov.incumb[[4]], lm.out.full.nocov.incumb.comp[[4]],
          #This next row includes the correct SEs:
          se = list(cl.se.full.nocov.incumb[[1]], cl.se.full.nocov.incumb.comp[[1]],
                    cl.se.full.nocov.incumb[[2]], cl.se.full.nocov.incumb.comp[[2]],
                    cl.se.full.nocov.incumb[[3]], cl.se.full.nocov.incumb.comp[[3]],
                    cl.se.full.nocov.incumb[[4]], cl.se.full.nocov.incumb.comp[[4]]),
          style = "ajps",
          label = paste("tab:lge-nocovariates-anconly"),
          column.separate = c(8),
          #covariate.labels = covars,
          title = "Heterogeneous effects of change in service delivery on change on local ANC vote: ANC municipalities only, No Covariates",
          model.numbers=F,
          dep.var.labels.include = FALSE,
          column.labels = c("Change in ANC Vote Share (Local Elections)"),
          #column.separate = c(1,1,1,1,1),
          omit.stat = c("f", "adj.rsq", "ser"),
          notes = "Standard errors clustered by municipality",
          float.env = "sidewaystable",
          font.size = "tiny",
          add.lines = list(c("Covariates?","N","N","N","N","N","N","N","N"),
                           c("$<$80th Percentile SD?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("ANC Control Only?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("Competitive only?","","Y","","Y","","Y","","Y")),
          out = paste(stgz,filename,sep="\\"))

###Table 19 in Appendix: LGE, Opp only, No Covariates
filename = paste("lgenocovariatesOPP.tex")
stargazer(lm.out.full.nocov.opp[[1]], lm.out.full.nocov.opp.comp[[1]],
          lm.out.full.nocov.opp[[2]], lm.out.full.nocov.opp.comp[[2]],
          lm.out.full.nocov.opp[[3]], lm.out.full.nocov.opp.comp[[3]],
          lm.out.full.nocov.opp[[4]], lm.out.full.nocov.opp.comp[[4]],
          #This next row includes the correct SEs:
          se = list(cl.se.full.nocov.opp[[1]], cl.se.full.nocov.opp.comp[[1]],
                    cl.se.full.nocov.opp[[2]], cl.se.full.nocov.opp.comp[[2]],
                    cl.se.full.nocov.opp[[3]], cl.se.full.nocov.opp.comp[[3]],
                    cl.se.full.nocov.opp[[4]], cl.se.full.nocov.opp.comp[[4]]),
          style = "ajps",
          label = paste("tab:lge-nocovariates-opponly"),
          column.separate = c(8),
          #covariate.labels = covars,
          title = "Heterogeneous effects of change in service delivery on change on local ANC vote: Opposition municipalities only, No Covariates",
          model.numbers=F,
          dep.var.labels.include = FALSE,
          column.labels = c("Change in ANC Vote Share (Local Elections)"),
          #column.separate = c(1,1,1,1,1),
          omit.stat = c("f", "adj.rsq", "ser"),
          notes = "Standard errors clustered by municipality",
          float.env = "sidewaystable",
          font.size = "tiny",
          add.lines = list(c("Covariates?","N","N","N","N","N","N","N","N"),
                           c("$<$80th Percentile SD?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("Opp Control Only?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("Competitive only?","","Y","","Y","","Y","","Y")),
          out = paste(stgz,filename,sep="\\"))

###Table 21 in Appendix: LGE, Competitive only
filename = paste("lgenocovariatesmain.tex")
stargazer(lm.out.full.nocov[[1]], lm.out.full.nocov.comp[[1]],
          lm.out.full.nocov[[2]], lm.out.full.nocov.comp[[2]],
          lm.out.full.nocov[[3]], lm.out.full.nocov.comp[[3]],
          lm.out.full.nocov[[4]], lm.out.full.nocov.comp[[4]],
          #This next row includes the correct SEs:
          se = list(cl.se.full.nocov[[1]], cl.se.full.nocov.comp[[1]],
                    cl.se.full.nocov[[2]], cl.se.full.nocov.comp[[2]],
                    cl.se.full.nocov[[3]], cl.se.full.nocov.comp[[3]],
                    cl.se.full.nocov[[4]], cl.se.full.nocov.comp[[4]]),
          style = "ajps",
          label = paste("tab:lge-no-covariates-main"),
          column.separate = c(8),
          #covariate.labels = covars,
          title = "Effect of change in service delivery on change on local ANC vote share, No Covariates",
          model.numbers=F,
          dep.var.labels.include = FALSE,
          column.labels = c("Change in ANC Vote Share (Local Elections)"),
          #column.separate = c(1,1,1,1,1),
          omit.stat = c("f", "adj.rsq", "ser"),
          notes = "Standard errors clustered by municipality",
          float.env = "sidewaystable",
          font.size = "tiny", 
          add.lines = list(c("Covariates?","N","N","N","N","N","N","N","N"),
                           c("$<$80th Percentile SD?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("Competitive only?","","Y","","Y","","Y","","Y")),   
          out = paste(stgz,filename,sep="\\"))

###Table 22 in Appendix: LGE, Covariates, ANC and Opp
filename = paste("lgecovariatesmain.tex")
stargazer(lm.out.full.cov[[1]], lm.out.full.cov.comp[[1]],
          lm.out.full.cov[[2]], lm.out.full.cov.comp[[2]],
          lm.out.full.cov[[3]], lm.out.full.cov.comp[[3]],
          lm.out.full.cov[[4]], lm.out.full.cov.comp[[4]],
          #This next row includes the correct SEs:
          se = list(cl.se.full.cov[[1]], cl.se.full.cov.comp[[1]],
                    cl.se.full.cov[[2]], cl.se.full.cov.comp[[2]],
                    cl.se.full.cov[[3]], cl.se.full.cov.comp[[3]],
                    cl.se.full.cov[[4]], cl.se.full.cov.comp[[4]]),
          style = "ajps",
          label = paste("tab:lge-covariates-main"),
          column.separate = c(8),
          #covariate.labels = covars,
          title = "Effect of change in service delivery on change on local ANC vote share",
          model.numbers=F,
          dep.var.labels.include = FALSE,
          column.labels = c("Change in ANC Vote Share (Local Elections)"),
          #column.separate = c(1,1,1,1,1),
          omit.stat = c("f", "adj.rsq", "ser"),
          notes = "Standard errors clustered by municipality",
          float.env = "sidewaystable",
          font.size = "tiny", 
          add.lines = list(c("$<$80th Percentile SD?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("Competitive only?","","Y","","Y","","Y","","Y")),   
          out = paste(stgz,filename,sep="\\"))

###Table 23 in Appendix: LGE, Covariates, ANC and Opp
filename = paste("lgecovariatesinteracted.tex")
stargazer(lm.out.full.interact.cov[[1]], lm.out.full.interact.cov.comp[[1]],
          lm.out.full.interact.cov[[2]], lm.out.full.interact.cov.comp[[2]],
          lm.out.full.interact.cov[[3]], lm.out.full.interact.cov.comp[[3]],
          lm.out.full.interact.cov[[4]], lm.out.full.interact.cov.comp[[4]],
          #This next row includes the correct SEs:
          se = list(cl.se.full.interact.cov[[1]], cl.se.full.interact.cov.comp[[1]],
                    cl.se.full.interact.cov[[2]], cl.se.full.interact.cov.comp[[2]],
                    cl.se.full.interact.cov[[3]], cl.se.full.interact.cov.comp[[3]],
                    cl.se.full.interact.cov[[4]], cl.se.full.interact.cov.comp[[4]]),
          style = "ajps",
          label = paste("tab:lge-covariates-interacted"),
          column.separate = c(8),
          title = "Heterogeneous effects of change in service delivery on change on local ANC vote share by local incumbent",
          #covariate.labels = covars,
          model.numbers=F,
          dep.var.labels.include = FALSE,
          column.labels = c("Change in ANC Vote Share (Local Elections)"),
          #column.separate = c(1,1,1,1,1),
          omit.stat = c("f", "adj.rsq", "ser"),
          notes = "Standard errors clustered by municipality",
          float.env = "sidewaystable",
          font.size = "tiny",
          add.lines = list(c("$<$80th Percentile SD?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("Competitive only?","","Y","","Y","","Y","","Y")),
          out = paste(stgz,filename,sep="\\"))

###Table 25 in Appendix: LGE, Excluding WC
filename = paste("lgecovariatescompnowc.tex")
stargazer(lm.out.full.cov.comp.nowc[[1]], lm.out.full.cov.nowc[[1]],
          lm.out.full.cov.comp.nowc[[2]], lm.out.full.cov.nowc[[2]],
          lm.out.full.cov.comp.nowc[[3]], lm.out.full.cov.nowc[[3]],
          lm.out.full.cov.comp.nowc[[4]], lm.out.full.cov.nowc[[4]],
          #This next row includes the correct SEs:
          se = list(cl.se.full.cov.comp.nowc[[1]], cl.se.full.cov.nowc[[1]],
                    cl.se.full.cov.comp.nowc[[2]], cl.se.full.cov.nowc[[2]],
                    cl.se.full.cov.comp.nowc[[3]], cl.se.full.cov.nowc[[3]],
                    cl.se.full.cov.comp.nowc[[4]], cl.se.full.cov.nowc[[4]]),
          style = "ajps",
          label = paste("tab:lge-covariates-nowc"),
          column.separate = c(8),
          #covariate.labels = covars,
          title = "Effects of change in service delivery on change on local ANC vote, Excluding Western Cape",
          model.numbers=F,
          dep.var.labels.include = FALSE,
          column.labels = c("Change in ANC Vote Share (Local Elections)"),
          #column.separate = c(1,1,1,1,1),
          omit.stat = c("f", "adj.rsq", "ser"),
          notes = "Standard errors clustered by municipality",
          float.env = "sidewaystable",
          font.size = "tiny",
          add.lines = list(c("Covariates?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("$<$80th Percentile SD?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("Competitive only?","Y", "", "Y","","Y", "", "Y", "")),
          out = paste(stgz,filename,sep="\\"))

###Table 27 in Appendix: LGE, Excluding KZN
filename = paste("lgecovariatescompnokzn.tex")
stargazer(lm.out.full.cov.comp.nokzn[[1]], lm.out.full.cov.nokzn[[1]],
          lm.out.full.cov.comp.nokzn[[2]], lm.out.full.cov.nokzn[[2]],
          lm.out.full.cov.comp.nokzn[[3]], lm.out.full.cov.nokzn[[3]],
          lm.out.full.cov.comp.nokzn[[4]], lm.out.full.cov.nokzn[[4]],
          #This next row includes the correct SEs:
          se = list(cl.se.full.cov.comp.nokzn[[1]], cl.se.full.cov.nokzn[[1]],
                    cl.se.full.cov.comp.nokzn[[2]], cl.se.full.cov.nokzn[[2]],
                    cl.se.full.cov.comp.nokzn[[3]], cl.se.full.cov.nokzn[[3]],
                    cl.se.full.cov.comp.nokzn[[4]], cl.se.full.cov.nokzn[[4]]),
          style = "ajps",
          label = paste("tab:lge-covariates-nokzn"),
          column.separate = c(8),
          #covariate.labels = covars,
          title = "Effects of change in service delivery on change on local ANC vote, Excluding KwaZulu-Natal",
          model.numbers=F,
          dep.var.labels.include = FALSE,
          column.labels = c("Change in ANC Vote Share (Local Elections)"),
          #column.separate = c(1,1,1,1,1),
          omit.stat = c("f", "adj.rsq", "ser"),
          notes = "Standard errors clustered by municipality",
          float.env = "sidewaystable",
          font.size = "tiny",
          add.lines = list(c("Covariates?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("$<$80th Percentile SD?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("Competitive only?","Y", "", "Y","","Y", "", "Y", "")),
          out = paste(stgz,filename,sep="\\"))

###Table 29 in Appendix: LGE, Excluding EC
filename = paste("lgecovariatescompnoec.tex")
stargazer(lm.out.full.cov.comp.noec[[1]], lm.out.full.cov.noec[[1]],
          lm.out.full.cov.comp.noec[[2]], lm.out.full.cov.noec[[2]],
          lm.out.full.cov.comp.noec[[3]], lm.out.full.cov.noec[[3]],
          lm.out.full.cov.comp.noec[[4]], lm.out.full.cov.noec[[4]],
          #This next row includes the correct SEs:
          se = list(cl.se.full.cov.comp.noec[[1]], cl.se.full.cov.noec[[1]],
                    cl.se.full.cov.comp.noec[[2]], cl.se.full.cov.noec[[2]],
                    cl.se.full.cov.comp.noec[[3]], cl.se.full.cov.noec[[3]],
                    cl.se.full.cov.comp.noec[[4]], cl.se.full.cov.noec[[4]]),
          style = "ajps",
          label = paste("tab:lge-covariates-noec"),
          column.separate = c(8),
          #covariate.labels = covars,
          title = "Effects of change in service delivery on change on local ANC vote, Excluding Easternn Cape",
          model.numbers=F,
          dep.var.labels.include = FALSE,
          column.labels = c("Change in ANC Vote Share (Local Elections)"),
          #column.separate = c(1,1,1,1,1),
          omit.stat = c("f", "adj.rsq", "ser"),
          notes = "Standard errors clustered by municipality",
          float.env = "sidewaystable",
          font.size = "tiny", 
          add.lines = list(c("Covariates?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("$<$80th Percentile SD?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("Competitive only?","Y", "", "Y","","Y", "", "Y", "")),
          out = paste(stgz,filename,sep="\\"))

###################################
######INCLUDING TOP QUINTILE#######
###################################
i = 1
out.mod.full.cov.incumb = list()
out.mod.full.cov.incumb.comp = list()
out.mod.full.cov.opp = list()
out.mod.full.cov.opp.comp = list()

lm.out.full.cov.incumb = list()
lm.out.full.cov.incumb.comp = list()
lm.out.full.cov.opp = list()
lm.out.full.cov.opp.comp = list()

for(iv in c("water_ch10", "flush_ch10", "refuse_ch10", "sd_means_ch10")){  
  dat.sub = dat[,c(iv,
                   "lgeanc_vs_ch10", "unempfrac_ch10", "log_income_ch10", "log_pop_ch10", "femalefrac_ch10", "blackfrac_ch10", 
                   "whitefrac_ch10", "indianfrac_ch10", "colouredfrac_ch10", "traddwell_ch10", "cat_b", "opp_control", "province2011",
                   "employment_broad2001","log_income2001","log_pop2001","sex2001","indian_frac2001","colored_frac2001","black_frac2001","white_frac2001",
                   "share_area_all_tbvc_ta","traddwell_proportion2001","refuse_perc2001","all_flush_perc2001","home_water_perc2001","sd_means2001",
                   "competitive")]
  dat.sub = na.omit(dat.sub)
  
  dat.sub.comp = dat.sub[which(dat.sub$competitive==1),]
  
  dat.sub.comp = na.omit(dat.sub.comp)
  
  iv.int = paste(iv,"*opp_control", sep="")
  
  if(iv=="water_ch10"){
    baseline = " + home_water_perc2001"
  } else if(iv=="flush_ch10"){
    baseline = " + all_flush_perc2001"
  } else if(iv=="refuse_ch10"){
    baseline = " + refuse_perc2001"
  } else {
    baseline = " + sd_means2001"
  }  
  
  ########REGRESSIONS:
  ##By incumbency status
  lm.out.full.cov.incumb[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub[which(dat.sub$opp_control==0),])
  lm.out.full.cov.incumb.comp[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub.comp[which(dat.sub.comp$opp_control==0),])
  
  lm.out.full.cov.opp[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub[which(dat.sub$opp_control==1),])
  lm.out.full.cov.opp.comp[[i]] = lm(paste("lgeanc_vs_ch10", " ~ ", iv.int, cov, baseline, sep=""), data=dat.sub.comp[which(dat.sub.comp$opp_control==1),])
  
  #######CLUSTER STANDARD ERRORS:
  out.mod.full.cov.incumb[[i]] = coeftest(lm.out.full.cov.incumb[[i]],vcov = vcovCluster(lm.out.full.cov.incumb[[i]], dat.sub[which(dat.sub$opp_control==0),]$cat_b))
  out.mod.full.cov.incumb.comp[[i]] = coeftest(lm.out.full.cov.incumb.comp[[i]],vcov = vcovCluster(lm.out.full.cov.incumb.comp[[i]], dat.sub.comp[which(dat.sub.comp$opp_control==0),]$cat_b))
  
  out.mod.full.cov.opp[[i]] = coeftest(lm.out.full.cov.opp[[i]],vcov = vcovCluster(lm.out.full.cov.opp[[i]], dat.sub[which(dat.sub$opp_control==1),]$cat_b))
  out.mod.full.cov.opp.comp[[i]] = coeftest(lm.out.full.cov.opp.comp[[i]],vcov = vcovCluster(lm.out.full.cov.opp.comp[[i]], dat.sub.comp[which(dat.sub.comp$opp_control==1),]$cat_b))
  
  i = i+1
}

###Generate appropriate SEs for the tables
cl.se.full.cov.incumb = lapply(out.mod.full.cov.incumb, `[`,,2)
cl.se.full.cov.incumb.comp = lapply(out.mod.full.cov.incumb.comp, `[`,,2)
cl.se.full.cov.opp = lapply(out.mod.full.cov.opp, `[`,,2)
cl.se.full.cov.opp.comp = lapply(out.mod.full.cov.opp.comp, `[`,,2)

###Appendix table 13: Main results including top quintile, ANC only, LGE
filename = paste("topqincludedlgecovariatesANC.tex")
stargazer(lm.out.full.cov.incumb[[1]], lm.out.full.cov.incumb.comp[[1]],
          lm.out.full.cov.incumb[[2]], lm.out.full.cov.incumb.comp[[2]],
          lm.out.full.cov.incumb[[3]], lm.out.full.cov.incumb.comp[[3]],
          lm.out.full.cov.incumb[[4]], lm.out.full.cov.incumb.comp[[4]],
          #This next row includes the correct SEs:
          se = list(cl.se.full.cov.incumb[[1]], cl.se.full.cov.incumb.comp[[1]],
                    cl.se.full.cov.incumb[[2]], cl.se.full.cov.incumb.comp[[2]],
                    cl.se.full.cov.incumb[[3]], cl.se.full.cov.incumb.comp[[3]],
                    cl.se.full.cov.incumb[[4]], cl.se.full.cov.incumb.comp[[4]]),
          style = "ajps",
          label = paste("tab:lge-covariates-anconly"),
          column.separate = c(8),
          #covariate.labels = covars,
          title = "Heterogeneous effects of change in service delivery on change on local ANC vote: ANC municipalities only",
          model.numbers=F,
          dep.var.labels.include = FALSE,
          column.labels = c("Change in ANC Vote Share (Local Elections)"),
          #column.separate = c(1,1,1,1,1),
          omit.stat = c("f", "adj.rsq", "ser"),
          notes = "Standard errors clustered by municipality",
          float.env = "sidewaystable",
          font.size = "tiny",
          add.lines = list(c("Covariates?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           
                           c("ANC Control Only?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("Competitive only?","","Y","","Y","","Y","","Y")),
          out = paste(stgz,filename,sep="\\"))

###Appendix table 14: Main results including top quintile, ANC only, LGE
filename = paste("topqincludedlgecovariatesOPP.tex")
stargazer(lm.out.full.cov.opp[[1]], lm.out.full.cov.opp.comp[[1]],
          lm.out.full.cov.opp[[2]], lm.out.full.cov.opp.comp[[2]],
          lm.out.full.cov.opp[[3]], lm.out.full.cov.opp.comp[[3]],
          lm.out.full.cov.opp[[4]], lm.out.full.cov.opp.comp[[4]],
          #This next row includes the correct SEs:
          se = list(cl.se.full.cov.opp[[1]], cl.se.full.cov.opp.comp[[1]],
                    cl.se.full.cov.opp[[2]], cl.se.full.cov.opp.comp[[2]],
                    cl.se.full.cov.opp[[3]], cl.se.full.cov.opp.comp[[3]],
                    cl.se.full.cov.opp[[4]], cl.se.full.cov.opp.comp[[4]]),
          style = "ajps",
          label = paste("tab:lge-covariates-opponly"),
          column.separate = c(8),
          #covariate.labels = covars,
          title = "Heterogeneous effects of change in service delivery on change on local ANC vote: Opposition municipalities only",
          model.numbers=F,
          dep.var.labels.include = FALSE,
          column.labels = c("Change in ANC Vote Share (Local Elections)"),
          #column.separate = c(1,1,1,1,1),
          omit.stat = c("f", "adj.rsq", "ser"),
          notes = "Standard errors clustered by municipality",
          float.env = "sidewaystable",
          font.size = "tiny",
          add.lines = list(c("Covariates?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           
                           c("Opp Control Only?","Y","Y","Y","Y","Y","Y","Y","Y"),
                           c("Competitive only?","","Y","","Y","","Y","","Y")),
          out = paste(stgz,filename,sep="\\"))
